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Rayleigh-Benard convection is studied and quantitative comparisons are made, where possible, 
between theory and experiment by performing numerical simulations of the Boussinesq equations 
for a variety of experimentally realistic situations. Rectangular and cylindrical geometries of varying 
aspect ratios for experimental boundary conditions, including fins and spatial ramps in plate sepa- 
ration, are examined with particular attention paid to the role of the mean flow. A small cylindrical 
convection layer bounded laterally either by a rigid wall, fin, or a ramp is investigated and our results 
suggest that the mean flow plays an important role in the observed wavenumber. Analytical results 
are developed quantifying the mean flow sources, generated by amplitude gradients, and its effect 
on the pattern wavenumber for a large-aspect-ratio cylinder with a ramped boundary. Numerical 
results are found to agree well with these analytical predictions. We gain further insight into the 
role of mean flow in pattern dynamics by employing a novel method of quenching the mean flow 
numerically. Simulations of a spiral defect chaos state where the mean flow is suddenly quenched is 
found to remove the time dependence, increase the wavenumber and make the pattern more angular 
in nature. 
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I. INTRODUCTION 

Rayleigh-Benard convection has played a crucial role 
in guiding both theory and experiment towards an un- 
derstanding of the emergence of complex dynamics from 
nonequilibrium systems H. However, an important miss- 
ing link has been the ability to make quantitative and 
reliable comparisons between theory and experiment. 

Nearly all previous three-dimensional convection cal- 
culations have been subject to a variety of limitations. 
Many simulations have been for small aspect ratios where 
the lateral boundaries dominate the dynamics, and as a 
result, complicate the analysis. When larger aspect ratios 
are considered, it is often with the assumption of peri- 
odic boundaries, which is convenient numerically yet does 
not correspond to any laboratory experiment. As a re- 
sult of algorithmic inefficiencies, or the lack of computer 
resources, simulations have frequently not been carried 
out for long times. This presents the difficulty in de- 
termining whether the observed behavior represents the 
asymptotic non-transient state, which is usually the state 
that is most easily understood theoretically. 

Fortunately, advances in parallel computers, numerical 
algorithms and data storage are such that direct numeri- 
cal simulations of the full three-dimensional time depen- 
dent equations are possible for experimentally realistic 



'Electronic address: 


Iipaul@caltech.edu 


caltech.edu/~stchaoE 





situations. We have performed simulations with exper- 
imentally correct boundary conditions, in geometries of 
varying shapes and aspect ratios over long enough times 
so as to allow a detailed quantitative comparison between 
theory and experiment. 

Alan Newell has made numerous important contri- 
butions to the discussion of pattern formation in non- 
equilibrium systems. In this paper, presented in this 
special issue in his honor, we give a survey of our recent 
results that touch on many of the issues he has raised, 
and in turn make use of some of the tools that he has 
helped develop to understand our simulations. 



II. SIMULATION OF REALISTIC 
GEOMETRIES 

We have performed full numerical simulations of the 
fluid and heat equations using a parallel spectral element 
algorithm (described in detail elsewhere |||). The veloc- 
ity u, temperature T, and pressure p, evolve according to 
the Boussinesq equations, 

a- 1 (d t + u«V) u = -S7p + RTz + V 2 u, (I) 

(dt + U'V^T = V 2 T, (2) 

V.m = 0, (3) 

where dt indicates time differentiation, z is a unit vec- 
tor in the vertical direction opposite of gravity, R is the 
Rayleigh number, and a is the Prandtl number. The 
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equations are nondimensionalized in the standard man- 
ner using the layer height h, the vertical diffusion time 
for heat t v = h 2 /n where k is the thermal diffusivity, 
and the temperature difference across the layer AT, as 
the length, time, and temperature scales, respectively. 

We have investigated a wide range of geometries in- 
cluding cylindrical and rectangular domains, which are 
the most common experimentally, in addition to ellipti- 
cal and annular domains. Rotation about the vertical 
axis of the convection layer for any of these situations is 
also possible but will not be presented here. All bounding 
surfaces are no-slip, u = 0, and the lower and upper sur- 
faces and are held at constant temperature, T(z = 0) = 1 
and T(z = 1) = 0. 

A variety of sidewall boundary conditions are shown 
in Fig. |]. Common thermal boundary conditions on the 
lateral sidewalls are insulating, n«VT = where h is a 
unit vector normal to the boundary at a given point, and 
conducting, T = 1 — z. In the future we will have the 
flexibility of imposing a more experimentally accurate 
thermal boundary condition by coupling the fluid to a 
lateral wall of finite thickness and known finite thermal 
conductivity that is bounded on the outside by a vacuum. 



Conducting Insulating 




FIG. 1: Four lateral sidewall boundary conditions utilized in 
the numerical simulations. The two thermal boundary condi- 
tions are conducting and insulating whereas the fin and ramp 
represent geometric conditions employed in experiments. 

In experiment, however, small sidewall thermal forc- 
ing can have a significant effect upon the resulting pat- 
terns and, as a result, finned boundaries have been em- 
ployed H [|, |5| . These are formed by inserting a very thin 
piece of paper or cardboard between the top and bottom 
plates near the sidewalls. This suppresses convection over 
the finned region (R ~ h 3 and the layer height has effec- 
tively been reduced) whereas in the bulk of the domain, 
i.e. the un-finned region, supercritical conditions prevail. 
This is accomplished numerically by extending a no-slip 
surface into the domain from the lateral sidewall. In all 
of our simulations we have chosen the vertical position of 
the fin to be z — 0.5 but this is not necessary. The result 
is that the supercritical portion of the convection layer 
is bounded by a subcritical region of the same fluid and 
hence with the same material properties. An additional 
effect is that the mean flow may extend into the finned 
region which presents an interesting scenario for explor- 



ing the effect of mean flows upon pattern dynamics that 
has been investigated both experimentally and theoret- 
ically by Pocheau and Daviaud || |j| and is discussed 
further below. 

The sidewalls can also have an orienting effect and 
ramped boundaries have been used as a "soft bound- 
ary" H in an effort to minimize this. By gradually de- 
creasing the plate separation as the lateral sidewall is ap- 
proached the convection layer eventually becomes critical 
and then increasingly subcritical. Using the spectral ele- 
ment algorithm we are able to investigate arbitrary ramp 
shapes: we have chosen to investigate the precise radial 
ramp utilized in recent experiments |?], ^| on a cylindrical 
convection layer. Again the mean flow is able to extend 
into the subcritical region. 

Perhaps the most common method employed exper- 
imentally to reduce the influence of sidewalls is to use 
a large aspect ratio T, where T — r/h in a cylindrical 
domain where r is the radius and r = L/h in a square 
domain where L is the length of side. Experiments can 
attain aspect ratios as large as ~ 500. However, the ma- 
jority of large aspect ratio experiments are for T < 100. 
We have performed numerical simulations using the spec- 
tral element algorithm for r ~ 60 as shown in Fig. |2[ 

The top panel in Fig. || illustrates the convection pat- 
tern present for the parameters of the classic paper 
where flow visualization was not possible. Although the 
simulation has only been performed for a short time 
tf ~ IOOtv, it appears that a slow process of domain 
coarsening jl0| is occurring. The bottom of Fig. || illus- 
trates the time dependent spatiotemporal chaotic state 
of spiral defect chaos These, and other, interest- 

ing large aspect ratio problems can now be addressed 
through the use of numerical simulation. 

Heuristically, using the spectral element algorithm on 
an IBM SP parallel supercomputer, it is our experience 
that it is practical to perform full numerical simulations 
for aspect ratios r ~ 30 for simulation times of tt ~ Th 
(36 hours on 64 processors), where 77, is the horizon- 
tal diffusion time for heat 77, = T 2 t v , and r ~ 60 for 
tf ~ 300t„ (36 hours on 256 processors) for e < 1, 
0-5 ^ c < 10, At w 0.01, and approximately cubic 
shaped spectral elements with an edge length of unity and 
11 th order polynomial expansions (where e = (R—R c )/R c 
and R c is the critical value of the Rayleigh number) . Of 
course for smaller domains the computational require- 
ments significantly decrease. 

A major benefit of numerical simulations is that a com- 
plete knowledge of the flow field is produced. For exam- 
ple, we have first used this to address a long standing 
open question concerning chaos in small cylindrical do- 
mains. The existence of a power-law behavior in the 
fall-off of the power spectral density derived from a time 
series of the Nusselt number was not understood || . The 
Nusselt number, N(t), is a global measurement of the 
temperature difference across the fluid layer. In cryo- 
genic experiments very precise measurements of N(t) are 
possible H however the flow field can not be 
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visualized easily. Subsequent room temperature experi- 
ments using compressed gasses allowed flow visualization 
at the expense of precise measurements of the Nusselt 
number || || |T|. 



to experiment because of the presence of the noise floor. 




FIG. 2: Numerical simulations of two large-aspect-ratio cylin- 
drical convection layers. The pattern is illustrated by con- 
tours of the thermal perturbation, dark regions represent 
cool descending fluid and light regions warm ascending fluid. 
Both simulations are initiated from random thermal pertur- 
bations and the lateral sidewalls are insulating. (Top) T = 57, 
a = 2.94, R = 2169.2 and t = 74t„. (Bottom) A spiral defect 
chaos state, r = 30, a = 1.0, R = 2950 and t = 254r„. 

By performing long-time simulations, on the order of 
many horizontal diffusion times, for the same parame- 
ters in cylindrical domains with a = 0.78 and for a range 
of e, with realistic boundary conditions, we had access to 
both precise measurements of the Nusselt number, Fig. [|, 
and flow visualization, Fig. ^, allowing us to resolve the 
issue Conducting sidewalls were used and all sim- 

ulations were initiated from small, ST ~ 0.01, random 
thermal perturbations. Flow visualization of the sim- 
ulations represented in Fig. [| display a rich variety of 
dynamics similar to what was observed in the room tem- 
perature experiments. Using simulation results, the par- 
ticular dynamical events responsible for the N(i) signa- 
ture were identified. The power-law behavior was found 
to be caused by the nucleation of dislocation pairs and 
roll pinch-off events. Additionally, the power spectral 
density was found to decay exponentially for large fre- 
quencies as expected for time-continuous deterministic 
dynamics. The large frequency regime was not accessible 
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FIG. 3: Plots of the dimensionless heat transport N(t) for 
reduced Rayleigh number e — 0.557,0.614,0.8,1.0,1.5, and 
3.0, labelled (i-vi) respectively (F = 4.72). For cases (i-v), 
At = 0.01, and for case (vi), At = 0.005 (At is the time 
step). 



III. ROLE OF MEAN FLOW 

The mean flow present in these flow fields, and in gen- 
eral for a < 1, plays an important role in theory |n], ||^] 
yet it is not possible to measure or visualize the mean 
flows in the current generation of experiments. In our 
simulations, however, we can quantify and visualize the 
mean flow. 

The mean flow field, U(x,y), is the horizontal veloc- 
ity integrated over the depth and originates from the 
Reynolds stress induced by pattern distortions. Recall- 
ing the fluid equations, Eqs. |l]) and (0), it is evident 
that the pressure is not an independent dynamic vari- 
able. The pressure is determined implicitly to enforce 
incomprcssibility, 



V 2 p = -<r -1 ^. (w.Vj 



Rd z T. 



(4) 



Focussing on the nonlinear Reynolds stress term and 
rewriting the pressure as p = p Q {x, y) +p(x, y, z) yields, 



p {x,y) 



dx'dy' ln(l/ \r-r'\ 



(5) 

In Eq. (0) the ln(l/|r — r'\) is not exact, in order to be 
more precise the finite system Green's function would 
be required. However, the long range behavior persists. 
This gives a contribution to the pressure that depends on 
distant parts of the convection pattern. The Poiscuillc- 
like flow driven by this pressure field subtracts from the 
Reynolds stress induced flow leading to a divergence free 
horizontal flow that can be described in terms of a verti- 
cal vorticity. 
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FIG. 4: Flow visualization showing the pattern (solid dark 
lines) and shaded contours of the vorticity potential, f, for 
e = 0.614 labelled ii) in Fig || (F = 4.72). Dark regions corre- 
sponding to negative vorticity generate clockwise mean flow 
and light regions to a positive vorticity generating a counter 
clockwise mean flow. The dark solid lines are zeros of the 
thermal perturbation at mid-depth illustrating the outline of 
the convection rolls. From top to bottom and left to right 
the panels are for t = 600, 605, 630, 650, 735, 785. The dislo- 
cations glide toward the right wall focus (shown here); during 
the next half period, the dislocations glide to the left wall 
focus. This left and right alternation continues for the entire 
simulation. 



The mean flow is important not because of its strength; 
under most conditions the mean flow is substantially 
smaller than the magnitude of the roll flow making it ex- 
tremely difficult to quantify experimentally. The mean 
flow is important because it is a nonlocal effect acting 
over large distances (many roll widths) and changes im- 
portant general predictions of the phase equation [ p0[ . 
The mean flow is driven by roll curvature, roll compres- 
sion and gradients in the convection amplitude. The re- 
sulting mean flow advects the pattern, giving an addi- 
tional slow time dependence. 

The mean flow present in the simulation flow fields, 
U s (x,y), is formed by calculating the depth averaged hor- 
izontal velocity, 



U a (x,y) 



u±(x, y, z)dz 



(6) 



where u± is the horizontal velocity field. Furthermore it 
will be convenient to work with the vorticity potential, 
C, defined as 



where oj z is the vertical vorticity and is the horizontal 
Laplacian. 

Six consecutive snapshots in time for the periodic dy- 
namics shown in Fig. || case ii) are illustrated in Fig. |^. 
One half period is displayed illustrating the nucleation of 
a dislocation pair and its subsequent annihilation in the 
opposing wall foci. The vorticity potential, £, is shown on 
a grey scale: dark regions represent negative vorticity and 
light regions represent positive vorticity which will gen- 
erate a clockwise and a counter clockwise rotating mean 
flow, respectively. The quadrupole spatial structure of 
Q in the first panel, i.e. four lobes of alternating posi- 
tive and negative vorticity with one lobe per quadrant, 
generates a roll compressing mean flow that pushes the 
system closer to a dislocation pair nucleation event. Dur- 
ing dislocation climb and glide the spatial structure of the 
vorticity potential is more complicated until the pan-am 
pattern is reestablished in final panel and a quadrupole 
structure of vorticity is again formed and the process re- 
peats. The dislocations alternate gliding left and right 
resulting is a slight rocking back and forth of the entire 
pattern with each half period which is visible in the dif- 
ferent pattern orientations in the first last panels. This 
alternation persists for the entire simulation. 

A numerical investigation of the importance of the 
mean flow for this small cylindrical domain was per- 
formed by implementing the ramped and finned bound- 
ary conditions. In all of these simulations the bulk region 
of constant R extended out to a radius tq = 4.72. In 
the finned case a fin at half height occupied the region 
4.72 < r < 7.66. In the ramped case a radial ramp in 
plate separation was given by, 



h(r) 



( l=HL n ) 

\ri— r J 



r <r 
r>r 



(8) 



vie 



(7) 



where r = 4.72, n = 10.0, and S = 0.15. 

The different mean wavenumber behavior (using the 
Fourier methods discussed in Jnj]) exhibited in these 
three different cases is shown in Fig. |[ As illustrated 
in Fig. |^ the behavior of the vorticity potential suggests 
an explanation. In the simulations with a rigid sidewall, 
not ramped or finned, the vorticity potential generates a 
mean flow that enhances roll compression, as described 
above. In the case of the finned and ramped boundaries 
the vorticity potential and the resulting mean flow are be- 
ing generated by gradients in the convection amplitude 
and are largely situated away from the bulk of the do- 
main. Furthermore, the mean flow generated is strongest 
in the subcritical finned or ramped region away from the 
convection rolls. This is demonstrated by comparing the 
average value of the mean flow over a fraction of the bulk 
of the domain, r < 1, where it was found that U s = 0.23, 
0.09, and 0.02 for the rigid, finned and ramped domains, 
respectively, and that the maximum flow field velocity is 
\u\ « 10. 

It is attractive to pursue the case of a radial ramp in 
plate separation because the variation in the convective 
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FIG. 5: Mean pattern wavenumber measurements for a cylin- 
drical convection layer, rn = 4.72, with rigid sidewalls (□), 
fin (o, ri — 7.66) and a spatial ramp in plate separation 
(o, n = 10.0, r c = 7.34, and S = 0.15). In all three cases 
the sidewalls are perfectly conducting and a = 0.78. For ref- 
erence, solid lines labelled E, N, and SV indicate the approx- 
imate location of the Eckhaus, Neutral and Skewed Varicose 
stability boundaries for an infinite layer straight parallel con- 
vection rolls. All patterns represented are time independent. 




amplitude caused by the ramp can be determined ana- 
lytically and the influence of a mean flow upon nearly 
straight rolls can be quantified f2l|| . Usually the mean 
flow can only be determined once the texture is known 
and it is hard to calculate because of defects acting as 
sources, in addition to the regions of smooth distortions. 

Near threshold an explicit expression for the mean 
flow, U, that advects the convection pattern is pfj 



U{x 7 y) = -7/cVx 



(k\Af) 



V ±Po {x,y) (9) 



where 7 is a coupling constant given by 7 = 0.42ct _1 (ct + 
0.34)(cr + 0.51) -1 , \A\ 2 is the convection amplitude nor- 
malized so that the convective heat flow per unit area 
relative to the conducted heat flow at R c is \A\ 2 R/R Cl 
p is the slowly varying pressure (see Eq. (||)) and Vj_ is 
the horizontal gradient operator. The vertical vorticity 
is then given by the vertical component of the curl of 
Eq. (B), 



Vj_ x U I = -7z«Vj_ x 



fcVi 



(k\A\ 2 )] . (10) 



Consider a cylindrical convection layer with a radial ramp 
in plate separation containing a field of x-rolls given by 
k = k Q x. The amplitude can be represented for large e Q , 
using an adiabatic approximation, as \A\ 2 = e(r)Jg for 
e > and \A\ 2 = for e(r) < as shown in Fig. M, mak- 
ing the amplitude a function of radius only \A\ 2 = f(r). 
Inserting \A\ 2 — f(r) into Eq. ( |Tc| ) yields, after some 



FIG. 6: Convection pattern and shaded contours of the vortic- 
ity potential, £, for a cylindrical convection layer, rn = 4.72, 
with rigid sidewalls, fin (ri = 7.66) and a spatial ramp in plate 
separation (ri = 10.0, r c = 7.34 and marked with a dashed 
line, and <5 = 0.15) shown top, middle and bottom, respec- 
tively. The convection pattern is illustrated by plotting zero 
contours of the thermal perturbation. In all three cases the 
sidewalls are perfectly conducting, a — 0.78 and R — 2804. 



manipulation, the following expression for the vertical 
vorticity, 



uJ. 



lk 



d 2 \A{ 
dr 2 



1 d\Af 
r dr 



sin 29. 



(11) 



The vorticity generated by the amplitude variation 
caused by the ramp is also shown in Fig. 0: there is 
a negative vorticity for r < r < r c and then a delta 
function spike of positive vorticity at r c . To correct for 
nonadiabaticity and to smooth |v4(r) | 2 near r c , the one- 
dimensional time independent amplitude equation [^2| is 
solved, 

= f (r)A + ^ 2 cos 2 9^~g o \A\ 2 A, (12) 

where£ 2 = 0.148, g a = (0.6995-0.0047cr- 1 +0.0083(j- 2 ) 
and e(r) is determined by 



e(r) 



e , r <r 

e (h 3 -h 3 c )/(l-hl), r>r 



(13) 
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Fig. |]cj. Again the simulation results compare well with 
theory even in the absence of adjustable parameters. 
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FIG. 7: A schematic illustrating the radial variation, for 
purely adiabatic conditions, of e (dashed line), |A| 2 (solid 
line), and uj z (solid line with arrow) for a cylindrical con- 
vection layer with a radial ramp in plate separation. Labelled 
ro and r c are the radial values where the ramp begins and 
where the ramp yields critical conditions, respectively. Note 
that for r > r c , \A\ 2 = in the adiabatic approximation. 



Simulation 




Theory 




(b) 



where h c = h(r c ). Equation (|12| ) is solved numerically 
using the boundary conditions d r A = at r — 0, and 
A = at r — rj . 

To compare these analytical results with simulation 
we have chosen to investigate a large-aspect-ratio cylin- 
der with a gradual radial ramp, defined by Eq. (|J), given 
by the parameters: r = 11.31, r*i = 20.0, S r = 0.036, 
and a = 0.87. For small e the amplitude A 2 (r) is unable 
to adiabatically follow the ramp, this nonadiabaticity re- 
sults in a deviation from e(r)/g as shown in Fig. |^a. 
However, as eo increases the amplitude A 2 (r) follows 
e(r)/g adiabatically almost over the entire ramp except 
for a small kink at r c . The structure of uj z depends upon 
this adiabaticity and is shown in Fig. |^b where we have 
used the solution to Eq. ( |l2| ) at 9 = 7r/4 in Eq. (11). 
This is not strictly correct since the non-adiabaticity of 
the amplitude is 9 dependent which will induce higher an- 
gular modes of the vorticity not given by Eq. (11). How- 
ever, the calculation should give a good approximation to 
the main sin 29 component of the vorticity. It is evident 
from Fig. ||b that the vertical vorticity, calculated from 
the simulation results as an angular average weighted by 
sin 29 has an octupole angular dependence (octupole in 
the sense of an inner and outer quadrupolc) and is well 
approximated by theory without any adjustable param- 
eters. 

The mean flow generated by these vorticity distribu- 
tions is determined by solving Eq. ( |l0| ) with the bound- 
ary condition C( r i) = 0. The vorticity potential is re- 
lated to the mean flow in polar coordinates by (U r , U$) = 
(r~ 1 dg(, —d r (). The vorticity potential is expanded ra- 
dially in second order Bessel functions while maintain- 
ing the sin 29 angular dependence. Of particular interest 
is the mean flow perpendicular to the convection rolls, 
U r {9 = 0) or equivalently U x (y = 0), which is shown in 
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Simulation 




Theory 


r ° (c) 



FIG. 8: Panel (a) shows the solution of Eq. @ plotted as 
A 2 (r), shown for comparison is e(r)/g - Panel (b) compares 
the vertical vorticity found analytically from Eq. (n2|) with an 
angular average, weighted by sin 26, of the vertical vorticity 
from simulation. Panel (c) compares the mean flow found 
analytically from Eq. (Q) with the mean flow from simulation. 
Parameters are ro = 11.31, r c — 13.20, ri = 20.0, S r = 0.036, 
a = 0.87 and e D = 0.025. 

To make the connection between mean flow and 
wavenumber quantitative it is noted that the wavenum- 
ber variation resulting from a mean flow across a field 
of x-rolls can be determined from the one-dimensional 
phase equation, 



Ud T (b = D\\d x 



(14) 



where the wavenumber is the gradient of the phase, k = 
8 x cb, L>j, = e o T~\ and t- 1 = lCl^a-fCSm)- 1 [@. 
Figure |9|a illustrates the wavenumber variation for a 
large-aspect-ratio simulation, k(r) for r < ro, and makes 
evident the roll compression, k(r — 0) > k(ro). Figure ^> 
compares the mean flow calculated from simulation to the 
predicted value of the mean flow required to produce the 
wavenumber variation shown in Fig. Oa. The agreement 
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is good and the discrepancy near tq, which is contained 
within one roll wavelength from where the ramp begins, 
is expected because the influence of the ramp was not 
included in Eq. This illustrates quantitatively that 

is in indeed the mean flow that compresses the rolls in 
the bulk of the domain. 



3.25 
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FIG. 9: Panel (a), the variation in the local wavenumber along 
the positive x-axis, or equivalently k(r) at 9 = 0. Panel (b), a 
comparison of the mean flow from simulation (solid line) with 
the predicted value (dashed line) calculated from Eq. ( p"l| ) 
using the wavenumber variation from panel (a). Simulation 
parameters, ro = 11.31, ri = 20, 5 r = 0.036, a = 0.87 and 
e = 0.171. 

Finally, to better understand the connection between 
mean flow and pattern dynamics, especially that of 
spatiotemporal chaotic states exhibiting both temporal 
chaos as well as spatial disorder, we apply a novel nu- 
merical procedure to eliminate mean flow from the fluid 
equation, Eq. (jl]), thereby evolving the dynamics of an 
artificial fluid with no explicit contributions from mean 
flow. In this way, we can then obtain quantitative com- 
parisons between the patterns generated by this artificial 
fluid with mean flow quenched and by the original fluid 
equation. 

We have applied this procedure to study spiral defect 
chaos (see bottom of Fig. ||) |0. Numerous attempts 
have been made to understand how a spiral defect chaos 
state is formed and how it is sustained. For example, 
experiments |^3|, have found that spirals transition to 
targets when the Prandtl number is increased. Owing 



to the fact that the magnitude of mean flow is inversely 
proportional to the Prandtl number, c.f. Eq. (^|), it was 
believed that spiral defect chaos is a low Prandtl number 
phenomenon for which mean flow is essential to their dy- 
namics. This is supported by studies of convection mod- 
els based on the generalized Swift-Hohenberg equation 
p5| , p6| |27j] , where spiral defect chaos is not observed un- 
less a term corresponding to mean flow is explicitly cou- 
pled to the equation. However, these observations are 
by themselves insufficient. For example, there are many 
other effects in the fluid equations that grow towards low 
Prandtl numbers, and there could be limitations in the 
Swift-Hohenberg modelling. We have applied our numer- 
ical procedure to this case to explicitly confirm the role 
of mean flow in the dynamics of spiral defect chaos. 




FIG. 10: Spiral defect chaos (left) and angular textures (right) 
obtained by quenching mean flow. The left panel is at t = 
152-Tj, and displays the pattern upon which the mean flow 
is quenched, the right panel is at t = 320r„. In both cases, 
R = 2950, a = 1.0 and the lateral sidewalls are insulating. We 
see that the spiral arms transition to angular textures when 
mean flow is quenched. Also, the quenched state is stationary. 

Recalling that we can approximate mean flow to be 
the depth-averaged horizontal velocity, c.f. Eq. (^J), we 
can first depth-average the horizontal components of the 
fluid equation, Eq. (|l|), to obtain a dynamical equation 
for the mean flow U s : 



a 1 9 t C/ s + (T 1 / dz(u*V)uj_ 
Jo 



-V± I dzp+V 2 L U s + [ dzd zz u±. (15) 
Jo Jo 

In this equation, the term — Vj_ Jq dzp can be absorbed 
into the nonlinear Reynolds stress term via Eq. (^|) and 
so will be ignored henceforth. The resulting equation is 
then a diffusion equation in U s with a source term F s = 
Jq 1 dz(u'V)u± — a J Q dzd zz u±. If this source term were 
not present, then U s , being the solution to a diffusion 
equation, evolves to zero with an effective diffusivity a, 
the Prandtl number. Thus, the role of F s is to act as a 
generating source for the mean flow U s . Subtracting it 
from the fluid equation, Eq. (|l|) , then results in the mean 
flow being eliminated. 
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In practice, we found that it is necessary to actually 
subtract F s multiplied by a constant to ensure that the 
magnitude of mean flow becomes zero. This can be un- 
derstood in terms of the necessity to correct for the fact 
that Eq. (||) is only an approximation to the flow field 
that advects the rolls given by 

U= [ dzg(z)u± (16) 
Jo 

where g(z) is a weighting function depending on the full 
nonlinear structure of the rolls. This is discussed further 
elsewhere psf . 

We have carried out this procedure by introducing the 
term F s to the right-hand-side of the fluid equation after 
a spiral defect chaotic state becomes fully developed, typ- 
ically after about one horizontal diffusion time starting 
from random thermal perturbations as initial condition. 
We see that the spirals immediately, on the order of a 
vertical diffusion time, "straighten out" to form angular 
chevron-like textures; see Fig. [u]. Unlike spiral defect 
chaos, these angular textures are stationary (with the ex- 
ception of the slow motion of defects such as the gliding 
of dislocation pairs). Thus, we have shown that when 
mean flow is quenched via the subtraction of the term 
F s from the fluid equation, spiral defect chaos ceases to 
exist. 

We have further quantified the differences between spi- 
ral defect chaos and the angular textures. We men- 
tion here briefly one of the results: by comparing the 
wavenumber distribution for both sets of states, we have 
observed that the mean wavenumber approaches the 
unique wavenumber possessed by axisymmetric patterns 
asymptotically far away from the center (29). (The ax- 
isymmetric pattern, by symmetry, does not have mean 
flow components.) We discuss this as well as other re- 
sults in a separate article j28|. 

IV. CONCLUSION 

Full numerical simulations of Rayleigh-Benard convec- 
tion in cylindrical and rectangular shaped domains for a 



range of aspect ratios, 5 < T < 60, with experimentally 
realistic boundary conditions, including rigid, finned and 
spatially ramped sidewalls, have been performed. These 
simulations provide us with a complete knowledge of the 
flow field allowing us to quantitatively address some in- 
teresting open questions. 

In this paper we have emphasized the exploration of 
the mean flow. The mean flow is important in a the- 
oretical understanding of the pattern dynamics, yet is 
very difficult to measure in experiment, making numeri- 
cal simulations attractive to close this gap. 

The mean flow is found to be important in small cylin- 
drical domains by investigating the result of imposing 
different sidewall boundary conditions. Analytical re- 
sults are developed for a large-aspect-ratio cylinder with 
a radial ramp in plate separation. Numerical results of 
the vertical vorticity and the mean flow agree with these 
predictions. Furthermore, the wavenumber behavior pre- 
dicted using the mean flow in a one-dimensional phase 
equation also agrees with the results of simulation. This 
allows extrapolation of the analysis to larger aspect ra- 
tios. 

Lastly we utilize the control and flexibility offered by 
numerical simulation to investigate a novel method of 
quenching numerically the mean flow. We apply this to a 
spiral defect chaos state and find that the time dependent 
pattern becomes time independent, angular in nature, 
and that the pattern wavenumber becomes larger. 

These quantitative comparisons illustrate the benefit of 
performing numerical simulations for realistic geometries 
and boundary conditions as a means to create quantita- 
tive links between experiment and theory. 
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